/*

   **********************************************************************
   *                                                                    *
   * This file contains the code for an image object for a content-     *
   * based image retrieval application that "summarizes" an image and   *
   * sends the resulting image "signature" to a server for              *
   * classification. The server returns the classification result in    *
   * the form of a Web page containing example images of the input      *
   * image's category.                                                  *
   *                                                                    *
   * The image class defined here is rather specific to this particular *
   * application and is not really designed for broad reuse (just yet). *
   *                                                                    *
   **********************************************************************
*/

/*
   **********************************************************************
   *                                                                    *
   * Gist   -- content-based image retrieval application for the gPhone *
   *                                                                    *
   * Copyright (C) 2011 Yong Zeng                                       *
   *                                                                    *
   * This file is part of Gist.                                       *
   *                                                                    *
   * Gist is free software; you can redistribute it and/or modify it  *
   * under the terms of the GNU General Public License as published by  *
   * the Free Software Foundation; either version 2 of the License, or  *
   * (at your option) any later version.                                *
   *                                                                    *
   * Gist is distributed in the hope that it will be useful, but      *
   * WITHOUT ANY WARRANTY; without even the implied warranty of         *
   * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU   *
   * General Public License for more details.                           *
   *                                                                    *
   * You should have received a copy of the GNU General Public License  *
   * along with Goggle; if not, write to the Free Software Foundation,  *
   * Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA.       *
   *                                                                    *
   **********************************************************************
*/


//----------------------- PACKAGE SPECIFICATION -------------------------

package vision.gistcomputing ;

//------------------------------ IMPORTS --------------------------------

// Android graphics support
import android.graphics.BitmapFactory ;
import android.graphics.Bitmap ;
import android.graphics.Color ;

// Android utilities
import android.util.Log ;

// Java I/O support
import java.io.FileOutputStream ;
import java.io.OutputStream ;
import java.io.InputStream ;

// Java utilities
import java.util.Arrays ;

//------------------------- CLASS DEFINITION ----------------------------

/**
   This class implements an image wherein the pixels are floating point
   values between 0 and 1 and are stored in a single-dimensional array
   organized in row-major order. No particular restrictions or
   assumptions are made regarding pixel formats. Thus, it is entirely up
   to clients whether they want to store RGB values or luminance values
   or whatever else they may need. In fact, this class is just a very
   thin wrapper around a float array (i.e., it is not really supposed to
   be some sort of high-level image API).
*/
class Image {

/// It is useful to keep track of the image's dimensions.
private int m_width, m_height ;

/// The image pixels are stored in an array of floats.
private float[] m_pixels ;

//----------------------------- ACCESSORS -------------------------------

/// Return the width of the image.
public int width()
{
   return m_width ;
}

/// Return the height of the image.
public int height()
{
   return m_height ;
}

/// Return the pixel value at the specified coordinates.
public float get(int x, int y)
{
   if (x < 0 || x >= m_width || y < 0 || y >= m_height)
      throw new IllegalArgumentException("bad coordinates: (" +
                                         x + ", " + y + ")") ;
   return m_pixels[y*m_width + x] ;
}

//-------------------------- INITIALIZATION -----------------------------

/// Returns a blank image of the specified size.
public Image(int width, int height)
{
   if (width <= 0 || height <= 0)
      throw new IllegalArgumentException("bad image size: " +
                                         width + "x" + height) ;

   m_width  = width ;
   m_height = height ;
   m_pixels = new float[width * height] ;
}

/// Copy constructor
public Image(Image I)
{
   this(I.m_width, I.m_height) ;
   System.arraycopy(I.m_pixels, 0, m_pixels, 0, I.m_pixels.length) ;
}

/// This function converts RGB color data to gray-scale and returns an
/// image with pixel values in the range 0-1.
public static Image gray_scale(int[] buffer, int W, int H)
{
   final int   N = W * H ;
   final float n = 1.0f/(3*255) ;

   Image I = new Image(W, H) ;
   for (int i = 0; i < N; ++i)
      I.m_pixels[i] = n * (  Color.red(buffer[i])
                           + Color.green(buffer[i])
                           + Color.blue(buffer[i])) ;
   return I ;
}

//------------------------- IMAGE PROCESSING ----------------------------

/**
   This function quadruples the size of the input image and fills in the
   extra pixels using linear interpolation (i.e., averaging) of
   neighbouring pixels from the source image.

   The process is as follows:

                                 a  A  b
                                 B  C  D
                                 c  E  d

   The above "diagram" shows a small patch of the newly created image.
   The lower-case letters indicate pixels from the original image whereas
   the upper-case letters are interpolated pixels. The interpolated
   pixels are computed as shown below:

                         A = (a + b) / 2
                         B = (a + c) / 2
                         C = (a + b + c + d) / 4
                         D = (b + d) / 2
                         E = (c + d) / 2

   The last row and column are processed specially using something
   half-way between a duplication and a black border.
*/
public static Image interpolate(Image I)
{
   if (I == null)
      return null ;

   Image II = new Image(I.m_width << 1, I.m_height << 1) ;

   final int sw =  I.m_width, sh =  I.m_height ;
   final int dw = II.m_width, dh = II.m_height ;

   float[] src = I.m_pixels ;
   float[] dst = II.m_pixels ;
   int s = 0,
      d = 0 ;

   // Interpolate all rows except the last one:
   for (int y = 0; y < sh - 1; ++y)
   {
      // Interpolate all columns on current row except last one:
      for (int x = 0; x < sw - 1; ++x)
      {
         dst[d]          = src[s] ; // top-left pixel
         dst[d + 1]      = 0.5f * (src[s] + src[s + 1]) ;  // top-right pixel
         dst[d + dw]     = 0.5f * (src[s] + src[s + sw]) ; // bottom-left pixel
         dst[d + dw + 1] = // bottom-right pixel
            0.25f * (src[s] + src[s + 1] + src[s + sw] + src[s + sw + 1]) ;

         ++s ;
         ++d ; ++d ;
      }

      // Interpolate last column of current row:
      dst[d]          = src[s] ;
      dst[d + 1]      = 0.75f *  src[s] ;
      dst[d + dw]     = 0.50f * (src[s + sw] + src[s]) ;
      dst[d + dw + 1] = 0.75f *  dst[d + dw] ;

      ++s ;
      d += dw + 2 ;
   }

   // Interpolate last row (except last pixel):
   for (int x = 0; x < sw - 1; ++x)
   {
      dst[d]          = src[s] ;
      dst[d + 1]      = 0.50f * (src[s] + src[s + 1]) ;
      dst[d + dw]     = 0.75f *  src[s] ;
      dst[d + dw + 1] = 0.75f *  dst[d + 1] ;

      ++s ;
      ++d ; ++d ;
   }

   // Finally, interpolate the bottom-right corner of the image:
   dst[d]          = src[s] ;
   dst[d + 1]      = src[s] * 0.75f ;
   dst[d + dw]     = src[s] * 0.75f ;
   dst[d + dw + 1] = src[s] * 0.50f ;

   return II ;
}

/// This function reduces the size of the input image by the given x and
/// y factors.
public static Image decimate(Image I, int X, int Y)
{
   if (I == null)
      return null ;

   // Don't go smaller than 1x1
   if (I.m_width <= 1 && I.m_height <= 1)
      return new Image(I) ;

   if (I.m_width == 1) // only need to thin out vertically
   {
      final int h = I.m_height/Y ;

      Image II = new Image(1, h) ;
      float[] src  =  I.m_pixels ;
      float[] dst  = II.m_pixels ;
      int s = 0, d = 0 ;

      for (int y = 0; y < h; ++y)
         dst[d++] = src[s++ + y * Y] ;

      return II ;
   }

   if (I.m_height == 1) // only need to thin out horizontally
   {
      final int w = I.m_width/X ;

      Image II = new Image(w, 1) ;
      float[] src  =  I.m_pixels ;
      float[] dst  = II.m_pixels ;
      int s = 0, d = 0 ;

      for (int x = 0; x < w; ++x)
         dst[d++] = src[s++ + x * X] ;

      return II ;
   }

   // Otherwise: decimate "regular" image (i.e., image whose dimensions
   // are greater than 1x1, 1xh or wx1).
   final int w = I.m_width/X ;
   final int h = I.m_height/Y ;
   final int skip = I.m_width % X + I.m_width * (Y - 1) ;

   Image II = new Image(w, h) ;
   float[] src  =  I.m_pixels ;
   float[] dst  = II.m_pixels ;
   int s = 0, d = 0 ;

   for (int y = 0; y < h; ++y)
   {
      for (int x = 0; x < w; ++x)
      {
         dst[d++] = src[s] ;
         s += X ;
      }
      s += skip ;
   }

   return II ;
}

/// This function reduces the size of the input image by a factor of
/// four.
public static Image decimate(Image I)
{
   return decimate(I, 2, 2) ;
}

/// This function applies a Gaussian filter with the specified standard
/// deviation to the input image.
public static Image blur(Image I, float sigma)
{
   return convolve(I, gaussian(sigma)) ;
}

/// This function blurs the input image using a Gaussian filter with a
/// standard deviation of 1.6.
public static Image blur(Image I)
{
   return blur(I, 1.6f) ;
}

/// This function removes salt-n-pepper noise from the given image and
/// returns the result as a new, denoised image. The denoising is
/// effected by first doubling the size of the image and filling in the
/// extra pixels using linear interpolation. Then, we apply a Gaussian
/// blur and, finally, decimate the image back to its original size.
public static Image denoise(Image I)
{
   /*
   Image interpolated = interpolate(I) ;
   Image blurred = blur(interpolated) ;
   Image decimated = decimate(blurred) ;
   //interpolated.save_jpeg("/data/data/cs546.group7/files/interpolated.jpg") ;
   //blurred.save_jpeg("/data/data/cs546.group7/files/blurred.jpg") ;
   //decimated.save_jpeg("/data/data/cs546.group7/files/decimated.jpg") ;
   return decimated ;
   //*/
   return decimate(blur(interpolate(I))) ;
}

//--------------- CONVOLUTION AND FILTERING OPERATIONS ------------------

/// This function convolves the input image with the supplied kernel. The
/// convolution operation is implemented as a separable filter wherein
/// the input image is first convolved with the kernel along the
/// x-direction. The result of this convolution is then convolved with
/// the kernel along the y-direction.
public static Image convolve(Image I, Image kernel)
{
   if (I == null)
      return null ;
   if (kernel == null) // return a copy of the input image
      return new Image(I) ;

   return yfilter(xfilter(I, kernel), kernel) ;
}

// This function convolves the input image with the supplied kernel along
// the x-direction.
private static Image xfilter(Image I, Image kernel)
{
   final int W  = I.m_width ;
   final int H  = I.m_height ;
   final int w  = kernel.m_width ;
   final int w2 = (w - 1)/2 ;

   // Flip the kernel (this makes the convolution a little more efficient)
   float[] K = new float[w] ;
   float kernel_sum = 0 ;
   for (int i = 0; i < w; ++i) {
      kernel_sum  += kernel.m_pixels[i] ;
      K[w - 1 - i] = kernel.m_pixels[i] ;
   }

   // Create a result image and setup "pointers" into the source and
   // destination pixel arrays.
   Image II = new Image(W, H) ;
   float[] src  =  I.m_pixels ;
   float[] dst  = II.m_pixels ;
   int s = 0, d = 0 ;

   // Do the convolution...
   if (W < w) // special case for small images (i.e., smaller than the kernel)
   {
      for (int y = 0; y < H; ++y)
         for (int x = 0; x < W; ++x)
         {
            float c = 0, n = 0 ;
            for (int k = 0; k < w; ++k)
               if (x + k - w2 >= 0 && x + k - w2 < W) {
                  c += src[s + k - w2] * K[k] ;
                  n += K[k] ;
               }
            dst[d++] = c * kernel_sum/n ;
            ++s ;
         }
   }
   else // convolution of regular sized images (i.e., larger than the kernel)
   {
      for (int y = 0; y < H; ++y)
      {
         // Left-most w2 pixels
         for (int x = w2; x > 0; --x)
         {
            float c = 0, n = 0 ;
            for (int k = x; k < w; ++k) {
               c += src[s + k - x] * K[k] ;
               n += K[k] ;
            }
            dst[d++] = c * kernel_sum/n ;
         }

         // The bulk of the pixels in the middle of the image
         for (int x = 0; x < W - w + 1; ++x)
         {
            float c = 0 ;
            for (int k = 0; k < w; ++k)
               c += src[s + k] * K[k] ;
            dst[d++] = c ;
            ++s ;
         }

         // Right-most w2 pixels
         for (int x = w - 1; x > w2; --x)
         {
            float c = 0, n = 0 ;
            for (int k = 0; k < x; ++k) {
               c += src[s + k] * K[k] ;
               n += K[k] ;
            }
            dst[d++] = c * kernel_sum/n ;
            ++s ;
         }

         s += w2 ;
      }
   }

   return II ;
}

// This function convolves the input image with the supplied kernel along
// the y-direction.
private static Image yfilter(Image I, Image kernel)
{
   final int W  = I.m_width ;
   final int H  = I.m_height ;
   final int h  = kernel.m_width ;
   final int h2 = (h - 1)/2 ;

   // Flip the kernel (this makes the convolution a little more efficient)
   float[] K = new float[h] ;
   float kernel_sum = 0 ;
   for (int i = 0; i < h; ++i) {
      kernel_sum  += kernel.m_pixels[i] ;
      K[h - 1 - i] = kernel.m_pixels[i] ;
   }

   // Create a result image and setup "pointers" into the source and
   // destination pixel arrays.
   Image II = new Image(W, H) ;
   float[] src  =  I.m_pixels ;
   float[] dst  = II.m_pixels ;
   int s = 0, d = 0 ;

   // Do the convolution...
   if (H < h) // special case for small images (i.e., smaller than the kernel)
   {
      for (int y = 0; y < H; ++y)
         for (int x = 0; x < W; ++x)
         {
            float c = 0, n = 0 ;
            for (int k = 0; k < h; ++k)
               if (y + k - h2 >= 0 && y + k - h2 < H) {
                  c += src[s + W * (k - h2)] * K[k] ;
                  n += K[k] ;
               }
            dst[d++] = c * kernel_sum/n ;
            ++s ;
         }
   }
   else // convolution of regular sized images (i.e., larger than the kernel)
   {
      // Since we're moving vertically along an image stored in row-major
      // order, consecutive elements will have their indices every W
      // spaces apart. Rather than compute these indices over and over
      // inside the convolution loops, we precompute them.
      int[] index = new int[h] ;
      for (int i = 0; i < h; ++i)
         index[i] = W * i ;

      // Top-most h2 pixels
      for (int y = h2; y > 0; --y)
      {
         for (int x = 0; x < W; ++x)
         {
            float c = 0, n = 0 ;
            for (int k = y; k < h; ++k) {
               c += src[s + index[k - y]] * K[k] ;
               n += K[k] ;
            }
            dst[d++] = c * kernel_sum/n ;
            ++s ;
         }
         s -= W ;
      }

      // The bulk of the pixels in the middle of the image
      for (int y = 0; y < H - h + 1; ++y)
      {
         for (int x = 0; x < W; ++x)
         {
            float c = 0 ;
            for (int k = 0; k < h; ++k)
               c += src[s + index[k]] * K[k] ;
            dst[d++] = c ;
            ++s ;
         }
      }

      // The bottom-most h2 pixels
      for (int y = h - 1; y > h2; --y)
      {
         for (int x = 0; x < W; ++x)
         {
            float c = 0, n = 0 ;
            for (int k = 0; k < y; ++k) {
               c += src[s + index[k]] * K[k] ;
               n += K[k] ;
            }
            dst[d++] = c * kernel_sum/n ;
            ++s ;
         }
      }
   }

   return II ;
}

//-------------------------- IMAGE GRADIENTS ----------------------------

/// An image's gradient is a vector field that shows, at each point, the
/// magnitude and direction of maximal brightness. Since these vectors
/// must be available at each point, the most convenient means of
/// representing the gradient vector field is to use an image that has
/// the same size as the input image. However, an image (as implemented
/// here) consists of only floating point scalar values and cannot
/// represent other quantities directly. We work around this limitation
/// by simply using two images: one for the magnitudes of the vectors at
/// each pixel and the other for the corresponding directions.
public static class Gradient {
   Image magnitude ;
   Image direction ;

   Gradient(Image mag, Image dir) {
      magnitude = mag ;
      direction = dir ;
   }
}

/// This function computes the gradient vector field of the input image
/// that shows, at each point, the magnitude and direction of maximal
/// brightness.
public static Gradient gradient(Image I)
{
   final int W = I.m_width ;
   final int H = I.m_height ;
   Image magnitude = new Image(W, H) ;
   Image direction = new Image(W, H) ;

   // Setup source and destination "pointers"
   float[] src  = I.m_pixels ;
   float[] mag  = magnitude.m_pixels ;
   float[] dir  = direction.m_pixels ;
   int s = 0, m = 0, d = 0 ;

   // First row of gradient vector field (GVF) is all zeros
   Arrays.fill(mag, m, m + W, 0) ;
   Arrays.fill(dir, d, d + W, 0) ;
   s += W ; m += W ; d += W ;

   // Compute the GVF for the inner rows
   for (int y = 1; y < H - 1; ++y)
   {
      // Left-most pixel is zero
      mag[m++] = dir[d++] = 0 ; ++s ;

      // Compute GVF for inner columns
      for (int x = 1; x < W - 1; ++x)
      {
         float dx = src[s + 1] - src[s - 1] ;
         float dy = src[s + W] - src[s - W] ;

         mag[m++] = (float) Math.sqrt(dx*dx + dy*dy) ;
         dir[d++] = (float) Math.atan2(dy, dx) ;
         ++s ;
      }

      // Right-most pixel is zero
      mag[m++] = dir[d++] = 0 ; ++s ;
   }

   // The last row of the GVF is all zeros
   Arrays.fill(mag, m, m + W, 0) ;
   Arrays.fill(dir, d, d + W, 0) ;

   return new Gradient(magnitude, direction) ;
}

//--------------------- IMAGE PROCESSING HELPERS ------------------------

/// This function generates a 1D Gaussian kernel with the specified
/// standard deviation and half-width.
public static Image gaussian(float sigma, int half_width)
{
   if (half_width <= 0)
      throw new IllegalArgumentException("bad half-width for Gaussian kernel: "
                                         + half_width) ;

   Image kernel = new Image(2 * half_width + 1, 1) ;
   float[] k = kernel.m_pixels ;

   k[half_width] = 1 ;
   final float _1_OVER_2_SIGMA_SQR = -0.5f/(sigma * sigma) ;
   for (int i = 1; i <= half_width; ++i)
      k[half_width + i] =
      k[half_width - i] = (float) Math.exp(i * i * _1_OVER_2_SIGMA_SQR) ;

   Utils.normalize(k) ;
   return kernel ;
}

/// This function generates a 1D Gaussian kernel with the specified
/// standard deviation and a half-width of 4 (so that the total size of
/// the resulting 2D Gaussian filter is 9x9).
public static Image gaussian(float sigma)
{
   return gaussian(sigma, 4) ;
}

//----------------------- MISCELLANEOUS HELPERS -------------------------

/// A convenience class to hold an image's width and height together in a
/// single structure.
public static class Dims {
   public int width, height ;
   public Dims(int w, int h) {width = w; height = h;}
}

//--------------------------- DEBUG SUPPORT -----------------------------

/// Dump a couple of pixels.
public void dump()
{
   StringBuffer s = new StringBuffer() ;
   final int n = Math.min(Math.min(m_width, m_height), 10) ;
   for (int i = 0; i < n; ++i)
      s.append(get(i, i) + " ") ;
   Log.e("Goggle", "Image.dump(): " + s.toString()) ;
}

/// Save the image to a JPEG file. We can then pull the image from the
/// phone and see whether or not it is alright or just garbled nonsense.
public void save_jpeg(String file_name)
{
   try
   {
      final int N = m_width * m_height ;
      int[] buffer = new int[N] ;
      for (int i = 0; i < N; ++i) {
         int color = (int)(m_pixels[i] * 255) ;
         buffer[i] = Color.rgb(color, color, color) ;
      }
      Bitmap B = Bitmap.createBitmap(buffer, m_width, m_height,
                                     Bitmap.Config.ARGB_8888) ;
      OutputStream os = new FileOutputStream(file_name) ;
      B.compress(Bitmap.CompressFormat.JPEG, 100, os) ;
      os.close() ;
   }
   catch (Exception e)
   {
      Log.e("Goggle", "unable to save image to " + file_name, e) ;
   }
}

//-----------------------------------------------------------------------

} 
